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The dynamical evolution of a Bose-Einstein condensate trapped in a one-dimensional lattice potential is in- 
vestigated theoretically in the framework of the Bose-Hubbard model. The emphasis is set on the far-from- 
equilibrium evolution in a case where the gas is strongly interacting. This is realized by an appropriate choice of 
the parameters in the Hamiltonian, and by starting with an initial state, where one lattice well contains a Bose- 
Einstein condensate while all other wells are empty. Oscillations of the condensate as well as non-condensate 
fractions of the gas between the different sites of the lattice are found to be damped as a consequence of the col- 
lisional interactions between the atoms. Functional integral techniques involving self-consistently determined 
mean fields as well as two-point correlation functions are used to derive the two-particle-irreducible (2PI) effec- 
tive action. The action is expanded in inverse powers of the number of field components A/", and the dynamic 
equations are derived from it to next-to-leading order in this expansion. This approach reaches considerably 
beyond the Hartree-Fock-Bogoliubov mean-field theory, and its results are compared to the exact quantum dy- 
namics obtained by A.M. Rey et al, Phys. Rev. A 69, 033610 (2004) for small atom numbers. 
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I. INTRODUCTION 

The dynamical evolution of ultracold atomic quantum gases 
driven far out of equilibrium are a subject of intensely grow- 
ing interest. Precision measurement techniques for many- 
body observables have been and are being developed with 
vigourous effort. This technology has triggered a strong de- 
mand for progress in theoretically describing non-equilibrium 
quantum many-body dynamics of strongly interacting sys- 
tems beyond mean-field approximations. Recent highlights of 
this development include, e.g., the variation and enhancement 
of the atomic interactions on the basis of Feshbach scatter- 
ing resonances (Hi EE Hi as well as the achievement of 
strongly correlated regimes within optical lattices |6i7]. With 
these techniques, e.g., a Bose-Einstein condensate can sud- 
denly be brought out of equilibrium, so that mean-field ap- 
proximations such as chosen in Gross-Pitaevskii and Hartree- 
Fock-Bogoliubov theory no longer suffice to describe its time- 
evolution. 

Mean-field approximations are valid as long as both, clas- 
sical statistical and quantum fluctuations are small. This is, in 
general, only the case in situations close to equilibrium. But 
even close to equilibrium, fluctuations can play an important 
role on long time scales. For example, the long-time evolu- 
tion towards the Bose-Einstein equilibrium state relies on the 
fluctuations appearing beyond mean-field order in a quantum 
many-body description. 

Moreover, kinetic descriptions beyond mean-field order 
which have been studied intensely in the past (in the context 
of cold atomic gases cf., e.g., Refs. H H [H El El El El 
n~5l El fl7ll and references therein), usually neglect the ini- 
tial dynamics directly after the change in the boundary condi- 
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tions which drive the system out of equilibrium. This short- 
coming is cured in dynamical approaches in which coupled 
equations of motion for the correlation functions are derived 
to describe the time evolution starting from a specific initial 
state. The built-up of correlations beyond mean-field order, 
in these equations, is usually taken into account by means of 
non-Markovian integrations over the evolution history of the 
correlation functions, see, e.g. ITlfl9ll2rll2lll22ll23ll. 

In this article we consider the time evolution of a Bose- 
Einstein condensate trapped in a one-dimensional lattice po- 
tential in the framework of the Bose-Hubbard model. The 
emphasis is set on the far-from-equilibrium evolution in a case 
where the gas is strongly interacting. This is realized by an ap- 
propriate choice of the parameters in the Hamiltonian. In the 
non-equilibrium initial state one lattice well contains a Bose- 
Einstein condensate while all other wells are empty. At NIST, 
such a system has recently been realized in experiment 1 24] . 

We find oscillations of the condensed as well as the non- 
condensed fractions of the gas between the different sites of 
the lattice. These oscillations are damped as a consequence of 
the collisional interactions between the atoms. Our theoretical 
approach is based on functional integral techniques involving 
self-consistently determined mean fields as well as two-point 
correlation functions. Specifically, the dynamic equations are 
derived from the two-particle-irreducible (2PI) effective ac- 
tion HEi I2H l27ll . Any approximations necessary in practice 
are made at the level of this action before the dynamic equa- 
tions are derived by functional differentiation. This ensures 
crucial symmetries like total particle number and energy to be 
automatically conserved. The 2PI effective action has been 
used to compare various perturbative approximations with the 
exact result for the dynamics of a strongly correlated one- 
dimensional lattice Bose gas |22]. The approach yields dy- 
namic equations which constitute a proper initial value prob- 
lem and therefore are in principle valid over all times, in- 
cluding initial and long-term evolution. In the Markovian 
limit, in which their validity is restriced to certain intermedi- 
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ate times, the well-known kinetic theory descriptions includ- 
ing the Gross-Pitaevskii and Bogoliubov-de Gennes equations 
as well as Landau and Beliaev quasiparticle damping can be 
rederived from these dynamic equations 

In this article we consider a non-perturbative approxima- 
tion which reaches substantially beyond the Hartree-Fock- 
Bogoliubov mean field theory. This non-perturbative ap- 
proach is based on a systematic expansion of the 2PI effective 
action in powers of the inverse number of field components Af 
OS]. We point out that this 1 /Af expansion of the 2PI ef- 
fective action goes substantially beyond and is not to be con- 
fused with the standard 1/Af expansion of the 1PI effective 
action. For example, it can be used to calculate critical expo- 
nents characterizing correlation functions in the vicinity of a 
phase transition, even in the limit of vanishing Af |29|. The 
2PI 1 /Af expansion has recently been used to describe the dy- 
namics of an ultracold atomic Bose gas far from equilibrium 
j23| • We compare our results to the exact quantum dynamics 
also obtained in Ref. 12211 for small atom numbers. 

The 2PI 1 /Af expansion to next-to-leading order yields dy- 
namic equations which contain direct scattering, memory and 
"off-shell" effects. It allows to describe far-from-equilibrium 
dynamics as well as the late-time approach to quantum ther- 
mal equilibrium. Recently, these methods have allowed im- 
portant progress in describing the dynamics of strongly in- 
teracting relativistic systems far from thermal equilibrium for 
bosonic 119tl30M3U l32l B3l 13411 as well as fermionic degrees 
of freedom lllljli. 

Our article is organized as follows: In Section HU we re- 
call the functional description of the quantum many-body dy- 
namics and present the set of coupled dynamic equations for 
the mean field and the two-point correlation functions. We 
discuss consequences of this distinction for the self energies 
which quantify the beyond-HFB contributions in the dynamic 
equations for the correlation functions. In Section [^] we 
present numerical results for the non-equilibrium dynamics 
in a one-dimensional lattice potential and compare these to 
the results from corresponding exact calculations. Our con- 
clusions are drawn in SectionllVl 



constant. Only local interactions, between atoms at one lattice 
site are taken into account. 



A. Generating functional for correlation functions 

The functional integral formulation of the quantum dynam- 
ics of a Bose-Einstein gas has been discussed extensively in 
the literature. For details we refer to our recent discussion l23ll 
of the non-perturbative effective action approach to far-from- 
equilibrium dynamics of an ultracold Bose gas. To simplify 
our notation we will, in the following, change to a represen- 
tation of the field ^ in terms of its real and imaginary parts, 
i.e., 



(2) 



such that \P = ($1 + 1^2) /V%- Furthermore, we include 
the lattice index n into the argument of the field i.e., 
x = (t,n) = (x Q ,n). 

The time evolution of quantum correlation functions can 
be derived from a generating functional. This involves fields 
$i obeying Bose commutation relations [&i(t, x), $j(t, y)\ = 
— &2,ij${x — y), where 02 denotes the Pauli 2-matrix, as well 
as the non-equilibrium normalized initial-state density matrix 
PD(t ): 



Z[J,K;p D ] = Tr 



pD(to)7~cexp\i / $i(x)Ji(x) 

1 Jx,C 

^ ^(a^foiO^-fo)}], (3) 



where f = L" (fao En 1 Summation over double indices, 
e.g. i = 1, 2, is implied. Tq denotes time ordering of the ex- 
ponential integral along the closed time path C which extends 
from the initial time xq = to the largest relevant time in the 
functional integral and back to xq = 0, i.e., our short hand 
notation means 



II. THE 2PI EFFECTIVE ACTION APPROACH TO THE 
DYNAMICS OF A LATTICE BOSE GAS 



L = /c <fa °? = (iT' fao+ X° <fco )?- <4) 



We study a non-relativistic system characterized by the 
Bose-Hubbard Lagrangian 



C(n,t) 



-j[*;(t)* n+1 (t) + * 
u 

~2 
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(i) 



Here, \& is a complex-valued, i.e. two-component, field de- 
fined in time t and at site n of the 1 -dimensional lattice. We 
use units where H = 1. J is the coupling which depends 
on the hopping probability and therefore on the lattice depth, 
e, denotes an external potential, and U a real-valued coupling 



Note that the closed-time contour can be seen as ensuring the 
normalization of the generating functional Z. The terms linear 
in the external fields Jj have been added to allow to generate 
correlation functions of order n by functional differentiation: 



5 n Z[J,K;p D ] 



i n SJ(xi) ■ ■■SJiXn) 



(5) 



J=K=0 



where the field indices have been suppressed. The fields ify 
are needed for the self-consistency condition of the two-point 
correlation functions as described below. 

The generating functional l|3} can be written, for Gaussian 
initial states, for which all correlation functions of order n > 3 
vanish at the initial time, as the functional integral (cf., e.g., 
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[37]) 



Here, 



Z[J,K} = 



Here, 



V<S> exp{i[s[$] + / $i(x)Ji(x) 

L ./x,C 



1 

+ - 



xy,C 



ix)K t] (x,y)^(y)]}. (6) 



S 2 S{<f>] 



8(f>i{x)8<t>j{y) 
2U 



8c{xo - 2/0 ) 



S[<&] 



~yjSc(x,y) (x)^ (x)$j (x)$, (x) 



+ <t> k (x)(j) k (x)8 lj /2\8 

It III 



(11) 



, (7) 



is the classical action depending on the fluctuating fields <£>;. 
The number of field components is A/" — 2, here. The delta 
function is defined as <5c(x, y) = 8 nm 8c(xo — j/q)> with x = 
(xq, n), y = (yo, m). The classical inverse propagator reads 

iD^(x,y) = -i5c(x,y)a 2: tjd yo 

- HiB(x,y)8e(x - yo)8ij, (8) 

where 



Hm{x,y) = —J(6 n +i, m + <5„, m+ i) +e n 8 n 



(9) 



is the one-body Hamiltonian. Note that the Gaussian initial 
conditions have been absorbed into the current terms by re- 
defining the fields J and K accordingly |37|. 



B. Effective action and dynamic equations 

In most situations the correlation functions up to order n — 
2 are of primary interest, i.e., the mean field 4>i{x) = ($j(x)) 
and the normal and anomalous one-body density matrices, 
which, in the basis relate to the diagonal (i — j) and 
off-diagonal (i = 3 — j) matrix elements of the two-point 
functions Gij(x,y) = (Tc<tj(x)<& j (y)) r , here generalized to 



unequal times xo 7^ t/n ll45ll . In Ref. 12 3H we have recalled 
the essential aspects of the 2PI (2-particle-irreducible) effec- 
tive action approach which allows in a transparent way to ob- 
tain equations of motion for these one- and two-point Green's 
functions within perturbative as well as non-perturbative ap- 
proximation schemes with respect to the interaction constant 
U. We summarize the main results in the following. 

Starting from the generating functional Z[J, K] = 
exp{iW[J, K}}, Eq. 0, one derives, by double Legendre 
transform of the generating functional W[J, K] of connected 
Green's functions, the 2PI effective action which is a func- 
tional of the one- and two-point functions 



is the classical inverse propagator, with x = (xo,n), y = 
(yo,m). r[</>, G] describes the quantum system completely, 
i.e., knowing it, one derives, by means of Hamilton's principle 
of stationary action, 



8V\M 

8(f>i(x) 



= 0, 



and 



ST[ct),G} 
SGij (x, y) 



= 0, (12) 



the exact many-body time evolution equations for <fi and G. 



As originally discussed in Refs. 125ll26l.l27ll T2[<j>, G] is rep- 
resented as the sum of all closed two-particle irreducible (2PI) 
diagrams involving only the field </>, the full propagator G and 
bare vertices oc U/Af B46I1 . To solve the dynamic equations 
approximations are made on the level of the effective action, 
by truncating the diagrammatic expansion of T2[4>,G\. An 
important advantage of this approach is that crucial symme- 
tries like total particle number, energy, momentum, angular 
momentum, etc., are automatically fulfilled irrespective of the 
particular truncation, cf., e.g. Ref. |23]. In the appendix, we 
provide an explicit expression for the total energy within the 
NLO l/M approximation employed in this work. 

The equations of motion resulting from Eqs. (I12> are most 
conveniently written in terms of the real valued statistical and 
spectral correlation functions, 



Fa(x,y) = 2<{*i( a; ).*j(v)})c. (13) 
Pij(x,y) = i([^i(x),^j(y)]) c , (14) 



respectively. These are related to G through 



Gtj(x,y) = Fij(x,y) - -py(x, y) sgn c (x - y ), (15) 



mc} = s[ 



-TrjlnG" 
2 1 



const. 



(10) 



where sgn c is the sign function which evaluates to 1 for 
x a > Vo an( l to — 1 otherwise. In this representation, the time 
ordering translates into time integration limits in the equations 



of motion: 

i&2,ijd XQ 



y.C 



2U \ 
-jjrF ij (x,x)J(f) j (x) 

(H 1B (x,y) + -r[<j) k (x)<j>k(x) 



Fkk{x,x)]5c(x,y)\4>i(x 



x 



(16) 



jf c [ - ia 2tik S c (x, z)d Z0 - Mik(x, z)] ( ^ j J jjj ) 



F kj(z,y) 

Pkj{z,y) 



dz ( 3{ X ' Z 

^ik 



Here, 



(17) 




FIG. 1: (Color online) Diagrammatic representation of the lead- 
ing order (LO) and next-to-leading order (NLO) contributions in 
the l/TV-expansion, to the 2PI part ^[0, G] of the 2PI effective 
action. The thick blue lines represent 2-point functions Gij(x,y), 
the red crosses field insertions <f>i(x), and the wiggly lines vertices 
USc(x, y). At each vertex, it is summed over double field indices i 
and integrated/summed over double time and space variables x. 



Mij(x,y) = 6 i3 H 1B (x,y) 



\6 k {x)4> k {x) + F kk (x,x)j8e(x,y) 



U 



(18) 



is the mean-field energy matrix which includes the in- 
dependent terms of the classical inverse propagator iG^\p 



Eq. M 11 . and the local part E^ (x) of the self energy 



2/ 



'SGij(x,y) 



(19) 



which, to derive the above equations, has been decom- 
posed into real and imaginary parts as Y<ij(x,y; <fi) = 

£g> (x)Sc (x,y) + E£ (x, y; 4>) ~ (x, y; 0) sgn c (x - 

2/0)0. 

We point out that, neglecting the right-hand sides of 
Eqs. ilbt , jilt , these equations constitute a set of time- 
dependent Hartree-Fock-Bogoliubov (HFB) equations for the 
mean field and the two-point functions, cf., e.g. Refs. I22LI23I1 . 
In this approximation, the dynamics of p decouples from 
that of </> and F. Neglecting also F, Eq. d!6i is the Gross- 
Pitaevskii equation. Given the exact self energy E, Eqs. 
and dl7> are the exact equations for the field cf> and the corre- 
lation functions F and p, respectively. Eq. dl7> is equivalent 
to the Schwinger-Dyson equation for the full Green's function 
G. 



C. Non-perturbative 1 /Af approximation 

To derive the quantum many-body time evolution, details 
about the self energy E are required, and these are, in general, 
only available to a certain approximation. In the following 
we will consider an expansion of T2, to next-to-leading order 
(NLO) in powers of the inverse number of field components 



Af d, IH IH which, in our special case, is TV = 2. In 
the context of a non-relativistic Bose gas, this approximation 
has been discussed in detail in Ref. 12311 . The contribution 
T2 [(f), G] to the effective action T, Eq. dlOi . then involves a 
leading (LO) and next-to-leading order (NLO) part which can 
be diagrammatically represented as shown in Fig. \l\ While 
the leading order contribution involves one diagram, in NLO a 
chain of bubble diagrams is resummed. All of these diagrams 
are proportional to the same power of 1/Af since each vertex 
scales with 1/Af, which is cancelled by the (blue) propagator 
loops which scale with Af since they involve a summation over 
the field indices from 1 to Af. Note that the Hartree-Fock- 
Bogoliubov (HFB) approximation is given by an action T 2 
which involves T^ and the first diagram of r^ LO in Fig.^ 
(b), cf., e.g., 12311 . We compare, in the next section, our results 
in NLO of the 1/ Af approximation with the HFB dynamics as 
well as with a perturbatively 1 47 ] reduced 1 /Af approximation, 
the second-order coupling expansion, which, besides the HFB 
diagrams, takes into account still the two diagrams of second- 
order in U, i.e. the 2nd and 5th diagrams in Fig.^(b) which 
contain two wiggly lines each. 

From T 2 [4>,G) = T^°[(f>,G} + T^i&G] we obtain, 



using Eq. dl9> , the self energies Eij(a:, y) 
(i/2)sgn c (x - 2/o)Ef,(x, y), with 



yF 



(x,y) 



^fj(x,y) 



2U_ 

'Jf 



A F (x,y) \A p (x,y) 
■\A p {x,y) A F (x,y) 



F i3 (x, y) 



(20) 



where A F<p {x,y) = If, p {x, y) + P F , p (x, y; I F . p ). Theresum- 
mation to NLO in 1/Af is expressed by the coupled integral 
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equations for I F . p |38]|48]: 



I F (x,y)\ = U 
I P (x,y) J J\f 



dz I p {x, z) 



F(x,y) 2 - \p{x,y) 2 
2F ij (x,y)p ij (x,y) 

F(z,y) 2 - y(z,y) 2 



2F ij (z,y)p ij (z,y) 
+ jf° dz ( \ 2F ij {z,y)p ij {z,y) 



at the initial time to = 0. Note that, using a general nPI ap- 
proach 1 37], also initial states with higher-order correlations 
can be taken into account. 

The total number of particles N to t,n(t) at lattice site n and 
time t is given as the sum of the number of condensed and 
excited atoms (x — (i, n)): 



(21) 



Here, F 2 — FijFij, etc. The functions P F , P , which contribute 
to p in the self energies ( I20> and vanish if (pi = 0, read 

2U 



N to t,n(t) = (¥(x)V(x)) = N c<n (t) + N c 
N cxc . n (t) = i [F u (x,x) - 1] , 



,(t), (25) 
(26) 

(27) 



P F (x,y;I F , p ) = -—r^H F (x,y) 



v.. 



where the constant —1 inside the parantheses takes into ac- 
count the zero-point fluctuations. Formally, this —1 originates 
from the spectral function p which, due to the Bose commu- 
tation relations, vanishes at equal times, except for 



+ / dz[H F (x,z)I p (z,y) + I F (x,z)H p (z,y)] 



dz [H p (x, z)I F (z, y) + I p (x, z)H F (z, y)] 



pi 2 (xo,n;x 0l m) = p 21 (x ,n; x ,m) = 5 n 



(28) 



(i 



By virtue of the 2PI effective action approach, total particle 
number as well as total energy are exactly conserved in all 
approximations obtained by truncating the diagrammatic ex- 
pansion of r 2 , cf.lXl 

We consider initial states where the normal and anoma- 
lous fluctuations are zero, i.e., all atoms are in the conden- 
sate: {&(x)*(y)) c = [F u (x,y) - 5 c (x,y)}/2 = 0, and 
+ 1 dv dwI F (x,v)H p (v,w)I p (w,y)j, (22) (%)£(y)) c = [F n (x, y) - F 22 (x, y)}/2 + iF 12 (x, y) = 0, 
' " for xq — yo — 0. Hence we choose 



Xo ryo 

dv / dw I p (x,v)H F (v,w)Ip(w,y) 
i) Jo 

x pv 

dv / dw Ip(x,v)H p (v 1 w)I F (w,y) 

>: JO 

Va ryo 



2U ( 

Pp(x,y;Ip, p ) = -jjr{H p (x,y) 



F n (0,n;0,m) = F 22 (0,n;0,m) = 5 nm /2, 
F 12 (0,n;0,m) = F 21 (0,n; 0,m) = 0. (29) 



dz [H p (x, z)I p {z, y) + I p (x, z)H p (z, y)] 

+ dv dw Ip(x,v)Hp(v,w)Ip(w 1 y)\, (23) 
Jya Jyo 

wherein the functions H F p are defined as 

H F (x,y) = -<j)i(x)F i: j(x,y)(t>j{y), 

H p (x,y) = -cj) l (x)p l: j(x,y)(l) j (y). (24) 

The technical procedure to solve the above dynamic equations 
in every time step requires the determination of the functions 
I(x, y), before the actual propagation of the respective corre- 
lation functions. 



III. FAR-FROM-EQUILIBRIUM TIME EVOLUTION OF A 
ONE-DIMENSIONAL LATTICE BOSE GAS 



A. Initial conditions 

In the following we present our results for the time evo- 
lution of a lattice Bose gas initially in a state far from equi- 
librium. The 2PI approach invoked in this work allows us to 
specify initial states which are Gaussian in the sense that only 
the mean-field and 2-point correlation functions are non-zero 



The initial values for pij are prescribed by Eq. (I28> . The con- 
densate fraction is chosen to be non-zero at a single site m 
only, 



<Aj(0,n) 



f 2N8 ll 5 n 



(30) 



where N is the total number of particles, and j = 1 implies 
that the initial mean field tp = {<f>i + i4> 2 )/\/2 is real. 



B. Numerical results 

The expansion in inverse powers of the number of field 
components Af allows for approximations which are inher- 
ently non-perturbative, i.e., it takes into account diagrams up 
to infinite order in the physical coupling constant U . Di- 
viding the equations of motion by the tunneling parameter 
J, neglecting quantum fluctuations! 49], an d taking into ac- 
count that the Fij(x, y) and 4>i{x)4>j (y) approximately scale 
with the total number of particles N, one finds that all terms 
in the NLO 1/J\f equations of motion dl6> . < l 1 71 scale with 
some power of NU /J. In the series representing the self en- 
ergy Ey , arbitrary high powers are resummed. The remain- 
ing terms which take into account the quantum fluctuations 
13^1 are reduced by powers of N since the spectral functions 
Pij(x,y) are of order one. Nevertheless, if NU/J > 1 one 
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FIG. 2: (Color online) Time evolution of an ultracold Bose gas in a 
1 -dimensional lattice with Nl = 2 lattice sites and periodic bound- 
ary conditions. The three columns show the same quantities, for a 
different total number of atoms N = 20, 40, and 80, respectively. 
U/ J is chosen such that the characteristic interaction parameter is 
NU I J = 4 throughout. Initially, all atoms are in a Bose-Einstein 
condensate at lattice site 1. The graphs in the first line show the total 
number of atoms iVtot,i at site 1, normalized to N. Due to number 
conservation, one has iVtot,2 = N — iVtot,i. The second and third 
lines show the condensate fractions N c<n /N at sites n = 1,2, and 
the fourth line adds these up, N c /N = ~}2 n N c ,n/N. The last line 
shows the occupation n q of the quasimomentum q = vk [in units 
of 1/Lattice spacing] 13 1 i . for v — 1, normalized to the total quasi- 
momentum occupation no + = N. Notice the different scale 
in the last line. Our results from the 1/N expansion to NLO are 
shown as a thick solid line. The (blue) dashed lines show the results 
from an exact calculation by A.M. Rey and coworkers, cf. Figs. 5, 
6 of Ref. 12211 . The (green) dash-dotted curves show the dynamics 
resulting from a 1/N expansion in perturbative 2nd-order coupling 
approximation (see text), the (red) dotted curves correspond to the 
Hartree-Fock-Bogoliubov approximation. With both sets of curves 
we reproduce the respective results in l22ll . 



expects that any perturbative or loop approximation of T2, 
neglecting terms with higher powers of NU/ J, fails. In or- 
der to probe the accuracy of the non-perturbative l/J\f ap- 
proximation we chose NU/ J to be larger than 1 and com- 
pare our results with the results of an exact numerical calcu- 
lation 1 22] for limited numbers of lattice sites and particles. 
Following Ref. 1 22] we have considered cases with a total 
number of N-^ = 2 and iV L = 3 lattice sites, and particle 
numbers between N = 8 and 80. Periodic boundary condi- 
tions were imposed. In Figs. [2]and[3]we present our results 
for the non-equilibrium quantum many-body evolution in the 
non-perturbative next-to-leading order (NLO) 1/N approxi- 
mation as well as the 'perturbative' Hartree-Fock-Bogoliubov 
(HFB) and second-order coupling approximations. In Fig.|2] 
we consider 2 lattice sites. The graphs in the first line show the 
total number of atoms iV t ot,i at site 1, normalized to N. Due 
to number conservation, one has iVtot,2 = N — iVtot,i- The 
second and third lines show the condensate fractions N cn /N 



at sites n = 1, 2, cf. Eq. (I26> . and the fourth line adds these 
up, N c /N — J2 n N c ,n/N. The last line shows the occupation 
n qu , of the quasimomentum q v = 2hv/Nl, v — 1 [q in units 
of 1 /Lattice spacing], normalized to the total quasimomentum 



occupation 



N. 



n, 



is defined as 



^ nm 

= 2W ^ { cos ^ n ~ TO )1 Fkk ( x > f ) + ^ ( x )fa (y) 

^ nm 

+ sm[q(n-m)] F 21 (x, y) - F 12 {x, y) 

+ M*)My) - Mx)Mv)] - 1}, (3D 

with x = {t,n), y — (t,m). In Fig. |5J we consider 3 lat- 
tice sites, and plot the fractions of condensed, excited, and 
total atom numbers for sites 1 and 2, as well as the total num- 
ber of condensed atoms iV c = iV Ci i + 2N C ^- We also plot 
the quasimomentum populations tiq and n 27T /3- Note that the 
respective populations of sites 2 and 3 are equal due to the 
periodic boundary conditions. For the quasimomentum popu- 
lations this means that the term proportional sm[q(n — m)] in 
Eq. (13 1 1 vanishes. 

The solid line shows, in all plots, the evolution according 
to the dynamic equations in NLO 1/M approximation. For 
comparison, we show, as a dashed line, the results of an exact 
calculation obtained by A.M. Rey and coworkers for a coher- 
ent initial state of N atoms at one site J22I1 . 

The dotted and dash-dotted lines show the corresponding 
dynamics in the HFB and second-order coupling approxima- 
tions, respectively, cf. Sect. Ill CI In the HFB approximation 
all non-local self-energy terms vanish in the equations of mo- 
tion takes \\6\ and ( I17K and fluctuations enter in the form of a 
local energy shift only. This approximation is fully classical, 
i.e., it does not take into account any quantum statistical fluc- 
tuations, cf., e.g., Ref. 13^1 for a detailed discussion. In the 
second-order coupling expansion of the NLO 1/N approxi- 
mation, non-local contributions from the self energy allow for 
damping. In Figs. |2 and |3 the dotted and dash-dotted lines 
precisely confirm the respective results of Ref. [22] for the 
HFB and 2 nd -order 1/N approximations. 

The time evolution shows different characteristic periods. 
At early times, the condensate oscillates coherently between 
the lattice wells. Only a small number of atoms is scattered 
from the condensate fraction into excited modes. This dy- 
namics is effectively described by a set of coupled non-linear 
Gross-Pitaevskii-like equations for the condensate mean field 
at the different lattice sites. These have the form of coupled 
single-particle Schrodinger equations, and the dynamics cor- 
responds to Rabi oscillations. The frequency of these oscil- 
lations is determined by the expansion parameter NU/ J. In 
the case of 2 lattice wells, the Bose-Hubbard lattice gas with 
periodic boundary conditions resembles a Josephson junc- 
tion with coupling energy Ej — 2JN, and charging energy 
E c = 2U, cf., e.g., I40ll4lll . For vanishing U, the system un- 
dergoes Rabi oscillations with frequency fl = 2Ej/N = 4 J, 
i.e., the period of these oscillations is T — (ir/2) J -1 . Note 
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FIG. 3: (Color online) Time evolution of an ultracold Bose gas in a 
1 -dimensional lattice with TVl = 3 lattice sites and periodic bound- 
ary conditions. Only site n — 1 is filled initially with N — 8 atoms. 
The interaction strength is U/J = 1/3, such that NU/J = 2.67. 
Besides the quantities already shown in Fig.|2| we show the depletion 
N e *c,n/N = (JVtot.n - N c ,n)/N at sites n = 1, 2. The populations 
of site n = 3 are the same as those at n — 2 due to the periodic 
boundary conditions. The curves indicate the same approximations 
as in Fig.|5| Notice the different scale in the third and fourth lines. 



also that the choice of parameters NU/J = 4, N 2 > 400 
sets the system into the Josephson regime N~ 2 (NU / J) <C 
(NU/J) <C N 2 . In this regime, the equilibrium state has a 
well defined relative phase between the sites, and small oscil- 
lations around this configuration can be described as a collec- 
tive excitation, the Josephson plasmon with plasma frequency 
wjp - y/Ej(E c + 4Ej/N 2 ) = 2J v /(JV!7/J)+4 EJ. We 
find, however, that the frequency of the large oscillations in 
Fig.|2]is closer to the Rabi frequency of an ideal gas than to the 
plasma frequency Wjp = Ay/2 J, for which Tjp ~ 0.357r/J. 
This is, of course, not contradictory since the derivation of the 
above expression for wjp requires the assumption of small os- 
cillations around equilibrium. 

Due to the interactions between the atoms, also higher- 
order classical statistical and quantum correlations build up. 
To leading order this means that atoms are exchanged be- 
tween the condensate and the non-condensate modes of the 
gas. These processes lead to a rapid destruction of the con- 
densate fraction. Hence, the coherence of the gas deteriorates 
which, on a somewhat longer timescale, leads to damping of 
the Rabi oscillations. 

All three approximations describe the dynamics well within 
the first period of coherent oscillations. This was expected 
since the system is initially mostly classical, with higher-order 
correlations and fluctuations playing a minor role. At larger 
times, all approximations fail to describe the dynamics ac- 
curately. While the purely classical HFB approxiation even 



qualitatively fails to show the correct damping behaviour, the 
higher-order approximations yield damping but quantitatively 
different results. Our comparisons for different atom num- 
bers N indicate that the accuracy of the NLO 1/Af approxi- 
mation improves with increasing N. All examples show that 
the condensate fraction is, at large times, underestimated in 
the second-order coupling expansion. Although it is overesti- 
mated in the NLO 1 /Af approximation, our results show that 
the non-perturbative 1/Af resummation is capable of taking 
into account high correlations which, at large times become 
relevant. 



IV. CONCLUSIONS AND OUTLOOK 

In this work we have studied the far-from-equilibrium dy- 
namics of an ultracold, strongly interacting one-dimensional 
lattice Bose gas on the basis of the 2PI effective action in next- 
to-leading-order (NLO) 1/Af approximation. The 2PI effec- 
tive action preserves, at any truncation, vital conservation laws 
as total particle number and energy, and allows to derive in 
a consistent way non-perturbative approximations which re- 
main applicable also for strong interactions. For weakly in- 
teracting systems close to equilibrium, the 2PI approach gives 
the well-known mean-field theories and their kinetic exten- 
sions including the dynamics and dissipation of small excita- 
tions. However, the 2PI approach goes substantially beyond 
the basic Hartree-Fock-Bogoliubov approximation and leads 
to a description of dynamical evolution far from equilibrium 
in the form of an initial-value problem. Moreover, it is not 
restricted to a small number of particles. 

The systems have been considered in the framework of the 
Bose-Hubbard model and consist of two and three lattice sites 
one of which is coherently populated initially while the others 
are empty. Tunneling through the barriers leads to Rabi-like 
oscillations between the lattice wells which are damped due 
to particle interactions. We have studied the time evolution of 
the condensate and excited fractions of the gas as well as of 
the quasimomentum distribution and compared these results 
with those from the exact, fully quantal calculation of the dy- 
namics from Ref. |22]. We find that in the non-perturbative 
1 /Af approximation the 2PI approach gives a reduced damp- 
ing due to higher-order correlations as seen in the exact so- 
lution. Quantitatively, the damping of the oscillations of the 
condensate fractions is underestimated in the NLO 1 /Af ap- 
proximation while it is overestimated for the excited fractions. 
The agreement with the exact results appears to improve when 
choosing a larger total number N of particles. Our findings 
are similar to those in Ref. 13411 . where the non-equilibrium 
dynamics of a single non-linear harmonic oscillator has been 
studied in the NLO as well as in a restricted NNLO 1/Af ap- 
proximation and compared with the exact evolution. As com- 
pared to these results, which were obtained in the symmetric 
phase where the mean field <f> vanishes, we do not find, in NLO 
1/Af, irregular behaviour at large times, like the revivals seen 
in the Hartree-Fock-Bogoliubov approximation. A systematic 
study of the /^-dependence of the accuracy of the 1/Af ap- 
proximation, extending the exact studies to a larger particle 
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number and lattice size, usin g, e .g., stochastic quantization 
techniques along the lines of 14211 . is the subject of ongoing 
work. We furthermore expect that the particular choice of an 
initial coherent-state population affects the ensuing damping 
behaviour within the Markov time. This requires a detailed 
study of the non-Markovian effects which shall be a further 
topic of a future publication. 
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APPENDIX A: NUMBER AND ENERGY CONSERVATION 

An important advantage of the 2PI effective action approach 
is that crucial symmetries like total particle number and en- 
ergy are automatically fulfilled irrespective of the particular 
truncation. As was shown in Ref. |23] number conservation 
is a consequence of the Noether theorem in conjunction with 
the invariance of the theory under transformations which are 
elements of the group 0(AT). Each 2PI diagram in the expan- 
sion of Y[<fi, G] separately carries this property such that any 
truncation of the series leads to number conservation. As a 
consequence, the total number 

N(t) = ^J[4>i(x)Mx)+Gii(x,x)] (Al) 

is conserved in time \23]. For the lattice gas, the spatial inte- 
gral means J g f(x) = /(*, n). 

Energy conservation follows from time translation invari- 
ance of r, cf., e.g., Ref. [43]. Consider the general transla- 
tions in continuous space and time which vanish at the bound- 
ary, x M — ► + e^(x), where £^(2;) is a time- and space- 
dependent infinitesimal 4-vector. The mean field and 2-point 
functions transform, under these translations, to leading order 
in e, as <f>i{x) — ► <f>i{x) + e" '{x)d*<j>i(x), and Gij(x,y) — > 
Gij(x,y) + e v (x)d^G i:j (x,y) + e v {y)d^Gij{x,y), respec- 
tively. Here, d£ = d/dx v , etc. One can show that under these 
transformations the variation of the 2PI effective action Y can 
be written as Y[<f), G] -> Y[(j>, G] + 5Y[cf>, G], with 



bYU.G 



J X 



Since, by virtue of the stationarity conditions H2i . the varia- 
tion SY vanishes for all solutions of the equations of motion 
for 4>i and Gij, an integration by parts shows that T^ v is the 
conserved Noether current for the time-space-translations: 



$Y[<f>,G} = - / e v {x)mT^{x)=Q. 



(A3) 

T^ u {x) is identified as the energy-momentum tensor, and the 
conservation law for total energy is expressed as d^T ti0 (x) = 
0ord t fd 3 xT 00 {t,x) = 0. 

To leading order in the 1/Af approximation we used the 
above space-time translations, and Eq. (IA2> . to derive the 
energy-momentum tensor. For the interaction terms which 
depend on the coupling U, however, as well as for the terms 
in NLO of the 1/JV approximation, it is more convenient to 
use a procedure known from field theory on curved space- 
time. For a space-time-dependent metric g pl/ (x), the energy- 
momentum tensor is defined as lE4ll 



T^(x) 



2 6Y[(f), G;g^] 



(A4) 



where y/—g denotes the square root of minus the determinant 
of g^. 

In the following we only quote the result for the total 
energy E(t) = J2 n Too(t,n) of the lattice Bose gas de- 
scribed by the Lagrangian Q. Using the above definitions 
and the useful relations Sg pL/ — ~g^ p g l/r7 Sg pa and 5y/—g = 
^J—gg^Sg^v/2, one obtains, in flat Minkowski space-time, 
with = diag{l, -1, -1, -1}: 



U 

2U r 
Jf 



ry 

~ (cj) 2 {x) + G H (x,x)) 2 ~ 2H(x,x) 



E i t ) = \ [ S c {x,y){H 1B {x 1 y)[(j) l {y)(t)i{x) +G ll {y,x)] 

" Jxy 

I{x, x) 

2i I H(x,z)I(x,z) 

I(x,z)H(z,u)I(u,x) j. (A5) 



Note that there is no integration over xq = t. Moreover, we 
use the definitions 



U 
Jf 



G 2 {x,y)~i / I(x,z)G 2 (z,y) 



H(x,y) = -c/) l (x)G lj (x,y)(j> : j(y), 



,(A6) 
(A7) 



)8?E v (a 



(A2) 



with G 2 = GijGij. The term proportional to [Gu(x, x)] 2 in 
Eq. (IA5> stems from the LO contribution in the expansion of 
r 2 in powers of 1/J\f, and the term I(x, x) as well as the last 
line from the NLO contributions, cf. Fig.Q 
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